/**
 * @author Gao Jian
 * @email gaoj@njupt.edu.cn
 * @create date 2021-06-03 10:23:40
 * @desc
 * 关于tiff掩膜数据的处理，包括连接破碎的小块成大块，删除块内空洞，较小的块
 */

#include "omp.h"
#include <iostream>
using namespace std;

// 取出等于指定值
void ain_regions_cpu(unsigned char *data, unsigned char value, unsigned sz) {
#pragma omp parallel for
  for (unsigned i = 0; i < sz; i++) {
    if (data[i] == value) {
      data[i] = 1;
    } else {
      data[i] = 0;
    }
  }
}

// 取出大于指定值
void bin_regions_cpu(unsigned char *data, unsigned char value, unsigned sz) {
#pragma omp parallel for
  for (unsigned i = 0; i < sz; i++) {
    if (data[i] >= value && data[i] != 255) {
      data[i] = 1;
    } else {
      data[i] = 0;
    }
  }
}

void dilateparts_cpu(unsigned char *output, unsigned char *input,
                     unsigned char value, unsigned int nw, unsigned int nh) {
#pragma omp parallel for
  for (unsigned i = 0; i < nw * nh; i++) {

    if (input[i] == value) {
      output[i] = value;
      continue;
    }

    unsigned iw = i % nw;
    unsigned ih = i / nw;

    unsigned char bbond[8] = {iw > 0,
                              iw<nw - 1, ih> 0,
                              ih<nh - 1, iw> 0 && ih > 0,
                              iw > 0 && ih < nh - 1,
                              iw < nw - 1 && ih > 0,
                              iw < nw - 1 && ih < nh - 1};
    unsigned int ids[8] = {i - 1,      i + 1,      i - nw,     i + nw,
                           i - 1 - nw, i - 1 + nw, i + 1 - nw, i + 1 + nw};

    for (unsigned char t = 0; t < 4; t++) {
      if (bbond[t] && input[ids[t]] == value) {
        output[i] = value;
        break;
      }
    }
  }
}

// 标记符号 = index+1
void init_regions_cpu(unsigned int *regions, unsigned char *data,
                      unsigned char value, unsigned sz) {
#pragma omp parallel for
  for (unsigned i = 0; i < sz; i++) {
    if (data[i] == value) {
      regions[i] = i + 1; // 标记符号 = index+1
    } else {
      regions[i] = 0; // 区域外标记为0
    }
  }
}

bool flip_regions_cpu(unsigned int *regions, unsigned char mode,
                      unsigned int nw, unsigned int nh) {
  unsigned nz[4] = {nh, nh, nw, nw};

  int ni[4] = {nw, nw, nh, nh};
  int ns[4] = {nw, nw, 1, 1};
  int oi[4] = {0, nw - 1, 0, nw * (nh - 1)};
  int dz[4] = {1, -1, nw, -nw};

  bool ba = false;
  if (mode < 4) {
    unsigned tnz = nz[mode], tni = ni[mode], toi = oi[mode], tdz = dz[mode],
             tns = ns[mode];

#pragma omp parallel for
    for (unsigned i = 0; i < tnz; i++) {
      for (unsigned id = toi + i * tns, n = 1; n < tni; id += tdz, n++)

        if (regions[id] != 0 && regions[id + tdz] != 0 &&
            regions[id] < regions[id + tdz]) {
          regions[id + tdz] = regions[id];
          if (!ba)
            ba = true;
        }
    }
  }
  return ba;
}

void area_regions_cpu(unsigned int *areas, unsigned int *regions, unsigned sz) {
#pragma omp parallel for
  for (unsigned i = 0; i < sz; i++) {
    areas[i] == 0;
  }

  for (unsigned i = 0; i < sz; i++) {
    if (regions[i] == 0)
      continue;
    areas[regions[i] - 1]++;
  }
}

void cleanparts_cpu(unsigned char *data, unsigned int *areas,
                    unsigned int *regions, unsigned char value_outer,
                    unsigned int threshold, unsigned sz) {
#pragma omp parallel for
  for (unsigned i = 0; i < sz; i++) {
    if (regions[i] == 0)
      continue;

    if (areas[regions[i] - 1] < threshold)
      data[i] = value_outer;
  }
}

int clean_parts_cpu(unsigned char *data, unsigned char val1, unsigned char val2,
                    unsigned int threshold, unsigned int nw, unsigned int nh,
                    const unsigned verbose = 0) {
  unsigned int sz = nw * nh;

  unsigned int *regions = new unsigned[sz];
  unsigned int *areas = new unsigned[sz];

  init_regions_cpu(regions, data, val1, sz);

  bool b0, b1, b2, b3;
  unsigned nit = 0;
  do {
    b0 = flip_regions_cpu(regions, 0, nw, nh);
    b2 = flip_regions_cpu(regions, 2, nw, nh);
    b1 = flip_regions_cpu(regions, 1, nw, nh);
    b3 = flip_regions_cpu(regions, 3, nw, nh);
    nit++;
    if (verbose)
      printf("+");
    if (nit % 50 == 0 && verbose)
      printf("  [%3d]\n  ", nit);
  } while (b0 || b1 || b2 || b3);
  if (verbose)
    printf("  [%3d]\n", nit);

  area_regions_cpu(areas, regions, sz);
  cleanparts_cpu(data, areas, regions, val2, threshold, sz);

  delete areas;
  delete regions;

  return 0;
}